Densities and energies of nuclei in dilute matter 
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We explore the ground-state properties of nuclear clusters embedded in a gas of nucleons with the 
help of Skyrme-Hartree-Fock microscopic calculations. Two alternative representations of clusters 
are introduced, namely coordinate-space and energy-space clusters. We parameterize their density 
profiles in spherical symmetry in terms of basic properties of the energy density functionals used 
and propose an analytical, Woods-Saxon density profile whose parameters depend, not only on the 
composition of the cluster, but also of the nucleon gas. We study the clusters' energies with the 
help of the local-density approximation, validated through our microscopic results. We find that the 
volume energies of coordinate- space clusters are determined by the saturation properties of matter, 
while the surface energies are strongly affected by the presence of the gas. We conclude that both 
the density profiles and the cluster energies are strongly affected by the gas and discuss implications 
for the nuclear EoS and related perspectives. Our study provides a simple, but microscopically 
motivated modeling of the energetics of clusterized matter at subsaturation densities, for direct use 
in consequential applications of astrophysical interest. 
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I. INTRODUCTION 

The understanding of the structure and evolution of 
neutron stars, and the cataclysmic events generating 
them, is strongly related to our knowledge of the micro- 
physics of dense matter. Additionally, the implementa- 
tion of the equation of state (EoS) in compact-star mod- 
eling presents a double challenge: First, knowledge of 
the structure and dynamics of hadronic matter over dif- 
ferent domains of temperature, density and composition 
is needed - including important phase transitions and in- 
teractions with leptons and photons; second, this knowl- 
edge must be weaved together and modeled in such a way, 
that numerical studies and simulations are feasible. Ter- 
restrial experiments on nuclei, including systems under 
extreme conditions, are performed and analysed world- 
wide in order to meet the former challenge. It remains a 
non-trivial task to model and tabulate the acquired in- 
formation into usable form for astrophysical applications. 

In terms of relevant degrees of freedom, matter at low 
temperature and at sub-saturation density is not highly 
exotic: it will primarily consist of nucleons, leptons, and 
photons and therefore belong in the domain of low-energy 
nuclear physics. Depending on the situation, baryonic 
matter may be almost isospin-symmetric or extremely 
neutron rich - as in the inner crust of a neutron star. 
A liquid-gas phase transition is expected to occur at 
sub-saturation densities, from non-uniform matter to un- 
bound nucleons. Of particular interest in the present 
work is the intermediate domain of clusterized matter, 
where nuclei - as it were, nuclear droplets - coexist with 
free nucleons, as well as leptons and photons, in ther- 
mal and chemical equilibrium. Again depending on the 
conditions, different nuclei may be found in this mixture: 
from very light ones, like deuterons and alpha particles, 
to iron-group nuclei, and much heavier or exotic species. 



A common approach to the composition of clusterized 
stellar matter, followed in the most widely used EoS's, 
namely those by J.M. Lattimer and F.D. Swesty [l| or 
H. Shen et al. 0, is the so-called single- nucleus approx- 
imation (SNA): Besides free nucleons, only one kind of 
light cluster (alpha particles) and one kind of heavy clus- 
ter are assumed to exist. The idea is to account in an 
average way for the properties of the statistical cluster 
distribution. The SNA may not affect very strongly ther- 
modynamical properties of matter at certain temperature 
and density domains of interest Q , but it may have con- 
sequences for dynamical processes dependent on reaction 
rates of specific nuclei [4[ and for the gas-liquid phase 
transition. Therefore, more modern approaches rely on 
an extended nuclear statistical equilibrium (NSE) con- 
cept, where the distribution of clusters over, in principle, 
all mass numbers is taken into account and obtained self- 
consistently in conditions of statistical equilibrium . 

Originally, the NSE was introduced to describe the re- 
action network taking place at the end of the evolution of 
massive stars in red supergiants. Being very diluted, nu- 
clei interact weakly and are almost not modified by the 
surrounding medium. These conditions naturally lead 
to the Saha equations. The NSE in the dense and hot 
matter in the core of supernovae was first applied in the 
EoS of Hillebrandt and Wolff @. In recent NSE imple- 
mentations the interactions between a cluster and 
the surrounding gas is treated in the so-called excludcd- 
volumc approach. The clusters and the gas of light par- 
ticles do not overlap in space and the clusters binding 
energy is kept as in the free limit. It is known, however, 
from quantal approaches, that the cluster properties are 



modified by the coexistence with a gas [13[. Such ap- 
proaches are presently limited to a system of light nuclei. 
It is clear that in-medium modifications of larger clusters 
would also change the composition of matter and the re- 
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latcd EoS. It is therefore important to take such effects 
into account. 

This work provides a microscopically based mod- 
elling of in-medium properties of nuclei, specifically 
their ground-state densities and energies. We focus on 
medium-mass and heavy clusters. Clusters are first stud- 
ied microscopically through Hartrcc-Fock calculations in 
Wigncr-Seitz cells. Clustering emerges from the quan- 
tal calculation since a dense cluster can be distinguished 
from a dilute phase. Two different types of clusters can 
however be defined: a representation in coordinate space, 
pertinent in excluded-volume approaches; and another 
representation in energy space, conceptually similar to 
existing quantal approaches for light clusters. Next, we 
study the energies of clusters embeded in a gas, with the 
help of the Local-Density Approximation (LDA), which 
we validate using our microscopic results. Thanks to 
our analytical model for density and the LDA, config- 
urations away from those accessible within Hartree-Fock 
or other variational calculations can be studied and the 
in-medium modified self-energies can be directly used in 
the EoS modelling at finite temperature. 

For the microscopic calculations we employ the SLy4 
Skyrme functional, but our in-medium densities and en- 
ergies explicitly depend on the isoscalar and isovector 
properties of the underlying interaction, and can thus 
readily be computed for different functionals. 

The outline of our paper is as follows. In Sec. HH ba- 
sic information is provided on our theoretical methods, 
namely the Skyrmc-Hartree-Fock model, the LDA and 
the Woods-Saxon functions. In Sec. IIII1 we analyse our 
microscopic results and define two types of clusters. In 
Sec. IIV1 we present and validate our analytical model 
for the density profiles of clusters in a gas. In Sec. |V] 
we study their energies within the LDA. We conclude in 
Sec. Ell 



THEORETICAL METHODS AND 
CONSIDERATIONS 

A. Microscopic calculations 



Our modclization of clusters embedded in a nucleon 
gas will initially be informed by microscopic Hartree-Fock 
calculations using Skyrme functionals. The Hartree-Fock 
equations are solved self-consistently in spherical sym- 
metry in a mesoscopic volume, which can be associated 
with a Wigner-Seitz cell. Its radius R ce u must be much 
larger than that of a typical nucleus. Here we will use 
Rceii = 35 fm. Implicit here is the assumption, that 
the clusters are sufficiently far from each other to not 
interact. The purpose of the mesoscopic volume in its 
present use is to allow for the formation of a uniform 
gas around the cluster. Its radius is not chosen via a 
condition of charge neutrality. It is simply an auxiliary 
quantity which will be eliminated from the final results. 

The equations are solved with mixed boundary con- 
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ditions: Dirichlct boundary conditions are imposed on 
odd single-particle states and von Neumann conditions 
on even states. 

In principle, Hartree-Fock-Bogolyubov calculations 
could be performed, in order to include pairing [l4j |. Su- 
perfluidity is an essential property of dilute star mat- 
ter [15l4l7| . but we expect the in-medium modifications 
due to pairing to be small in comparison with the in- 
medium modifications due to the mean-field. As we will 
show, these latter can be recasted into in-medium mod- 
ified bulk and surface terms of the cluster energy func- 
tional. To correctly constrain such terms from the micro- 
scopic calculations and disentangle them from possible 
pairing terms, we therefore prefer to neglect pairing and 
stick to a Hartree-Fock modeling. 

The ener gy d ensity for a traditional Skyrme functional 
is given by [1S|] 



H(r) 



2m r 



2m„ 

+C p 2 + C 3P a+2 + C cSP t 
+D pl + D zP a p\ + D cB p 3 T 3 
+Ci 2 pV 2 p + L»12p 3 V 2 p 3 

-^°[3pV- J + p 3 V- J 3 ], 



(1) 



where 



Co — gio i C3 — TjLf 3 

C ctf = ±[3^ + t 2 (4x 2 + 5)} 

D = -\t {2x + 1) i D 3 = -±t 3 (2x 3 + 1) 

D cS = ±[-h{2 Xl + 1 ) + t 2 (2x 2 + 1)] 

Ci 2 = ^[-9^+^(5 + 4x2)] 

D 12 = ±[3t 1 (±+x 1 )+t 2 (i+x 2 )] 

and ti, xi are the usual Skyrme parameters. In the above, 

P = Pp + Pn ! P3 = Pn - Pp (2) 

are the isoscalar and isovector densities, while 

T = T p +T„ J T 3 — T n ~ T p (3) 

give the isoscalar and isovector expectation-value densi- 
ties of the momentum-squared operator and 



J — Jn 



Jn 



J3 — J n Jp 



(4) 



are the current densities. The densities are functions of 
the radial variable r and the subscript p or n specifies 
the densities of protons or neutrons, respectively. 

The SLy4 Skyrme functional [H, EH has been chosen 
at present for the microscopic calculations. This func- 
tional is widely used to calculate ground and excited 
states of nuclei and it is particularly interesting in studies 
of dense matter, because is has been adjusted to a realis- 
tic EoS of symmetric and neutron matter, obtained vari- 
ationally using the UVI4 and UVII potentials (2ll 122]. 
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We stress, however, that our modeling will be general 
in scope, so that it can be easily applicable to different 
functionals and interactions. 

In practice we will study isotopic chains ranging from 
Z = 20 (Ca) to Z = 126 (Pb). For each Z we vary 
the neutron number from somewhat below Z to approxi- 
mately 3000. We thus obtain a variety of solutions, from 
stable nuclei to very neutron-rich clusters embedded in a 
neutron gas of density up to, approximately, 0.02 frn -3 . 
We are therefore able to study microscopically neutron- 
rich clusterized matter, before the onset of deformation 
and pasta phases. The energies and density profiles ob- 
tained from the microscopic calculations will be used to 
validate our analytical models. 

The Coulomb interaction among protons in stellar 
matter is screened by the neutralizing electron back- 
ground. The associated modification to the cluster 
Coulomb energies is an important in-medium effect which 
is already accounted for in all EoS-NSE based mod- 
els [1, [l0[ . A possible extra effect of Coulomb screen- 
ing on the nuclear part of the energy functional through 
a modification of the density profile is expected to be 
small in the very neutron rich nuclei that constitute our 
microscopic sample [23}. 



B. Energy in the local-density approximation 

Within the Local-Density Approximation (LDA), the 
local energy density is determined by the proton and neu- 
tron densities and their gradients, while current densities 
are ignored. The LDA is closely related to the Thomas- 
Fermi approach whereby the ground-state properties can 
be determined by minimizing the energy with respect to 
variations of the density. For our purposes, using the 
LDA will transpose the problem of determining the clus- 
ter energy in-medium correction towards the determina- 
tion of the nuclear density profile in the Wigner-Seitz 
cell. 

Within the LDA, the kinetic energy density is also 
written in terms of the local density, being replaced by 
its value in nuclear matter of the corresponding density. 
In spherical symmetry we thus write, 



(5) 



where kp p ^ n = ("Sit 2 pp^n) 1 ^ is the Fermi momentum of 
protons and neutrons in infinite matter of density p p , n - 
With these considerations, the energy density in the 
LDA H]) is given by 

n LDA (r) = 



2 m, p 2m ri 



+C oP 2 + C 3 p a+2 + C cSP t 
+D pj + D 3 p a P 2 + D cSP3 t 3 
+C 12P V 2 p + D 12 p 3 \7 2 p 3 , (6) 



where 



(3- 2 ) 2/3 P^ 



(7) 



and r, t 3 are given by Eq. ([3]) . The total energy in a 
spherical cell of radius R ce ii will be given simply by 



Eld, 



An 



U L DA(r)r 2 dr. 



(8) 



All that is needed for computing the energy of the system 
is then its density profile - both the isoscalar and isovector 
components. 



C. Nuclear density profiles 

Density profiles of medium mass and heavy nuclei are 
known to be well described by Woods-Saxon profiles, 



p WS (r) 



Po 



l + exp[(r - R ws )/a] 



(9) 



where R w is the radius and a is the diffuseness param- 
eter of the Woods-Saxon profile. Other analytical forms 
have been considered in the literature [HI, [25j] and re- 
sults marginally depend on the chosen form. We note 
that Hartree-Fock calculations yield density profiles with 
small ripples in the bulk, but those are to a large extent 
artefacts of the mean-field approach and are expected to 
be washed out by correlations, with the obvious excep- 
tion of exotic bubble-shaped nuclei [26j] . 

A straighforward generalization of the above expres- 
sions for nuclei in the presence of a homogeneous gas 
with density p gas in spherical symmetry is given by: 



Po 



1 + cxp [(r-R ws )/a) 

j Pgas 

l + cxp[-(r~ R ws )/a 

PO — Pgas 



1 + cxp [{r-R ws )/a\ 
It is straightforward to show that 

PO — Pgas 



Pgas 



Ap'{R ws ) 



and 



€i{R ws ) 



PO + Pgas 



(10) 



(11) 



(12) 



Analogous expressions may be introduced separately for 
the proton and neutron densities. 

The density profiles of ordinary nuclei existing on earth 
or produced in nuclear facilities have been analyzed be- 
fore in semi-classical Thomas-Fermi models in terms of 
Woods-Saxon functions (see for instance Ref. [13]), cor- 
responding in our case to setting the external density 
Pg as to zero. Application of the Extended Thomas- 
Fermi model for dilute nuclei in neutron stars has also 
been considered on the basis of Woods-Saxon functions 
in Ref. [28| . In the case where only one kind of particles 
is dripping out of the cluster, the neutrons for instance, 



4 



the external proton density p ga s,p is zero. In this work 
we mostly analyze neutron rich nuclei where p ga s,p = 0. 
However, our approach can also be applied in cases where 
there are protons in the gas. 



III. ANALYSIS OF THE DENSITY PROFILE 

In this Section we analyze the density profiles obtained 
microscopically and we define two kinds of representa- 
tions for the in-medium clusters. Simple relations be- 
tween these two representations will also be derived. At 
first we fit the microscopic Hartree-Fock profiles with 
Woods-Saxon functions, in order to determine useful pa- 
rameters, e.g., the bulk density and asymmetry. 

The microscopic density profile of nucleonic matter 
Pceiiir) is represented on the top part of Fig. [T] in the 
specific case of 1950 Sn (triangles). On the bottom part 
of the figure we show also the profile of the number of 
particles 4:irp ce u (r)r 2 dr from which is obtained the total 
number of particles in the cell, 



-4, 



47T 



r 2 p C eii(r) dr. 



(13) 



It is clearly visible from the top part of Fig. [T] (left and 
right) that there is a dense nuclear cluster in the center 
of the cell and a dilute neutron gas in its outer part. 
However, the central cluster contains only a fraction of 
the number of particles as it is shown on the bottom part 
of Fig. []] and in this case, most of the particles are in the 
non-central region. 

There is an ambiguity on how to define dilute clusters, 
related to the overlap of the cluster with the gas: does the 
external gas penetrate inside the cluster, or is it excluded 
from the cluster? These two representations of the same 
system are depicted in Fig. [T] where the two panels on 
the left show how to build up the microscopic cell den- 
sity pceiiif) from two excluded sub-densities, while the 
two panels on the right show how the penetration of the 
external gas inside the cluster can be taken into account. 
We name these two representations the coordinate-space 
clustering (left of Fig. [1} and the energy-space clustering 
(right of Fig.H}, respectively. The coordinate-space clus- 
tering is the representation of clusters which naturally 
emerges if the density is taken as the relevant degree of 
freedom, as in the density functional theory. In this pic- 
ture, clusters are recognized as density fluctuations on 
top of a homogeneous medium, and occupy a volume in 
space. 

In the following, we aim at providing a global model 
for the density profile, but first, we shall discuss the two 
types of clustering. 



A. Coordinate-space clustering 

Clusters defined in coordinate space occupy a volume 
in space and are surrounded by a gas of light particles. 



In the following, these clusters will be referred to as the 
r-clusters and the gas as the r-gas. The frontier between 
the r-cluster and the r-gas is located in space at the place 
where the nucleonic density changes drastically. Here- 
after, it will be referred to as the surface of the r-cluster. 
Since the r-cluster has a surface diffuscness, the frontier 
cannot be sharp. At the surface of the r-cluster, there 
is therefore a narrow region in space where the r-cluster 
and the r-gas overlap, with a typical size of 1-2 fm, as 
illustrated on the left part of Fig. [TJ 

The coordinate-space clusters have been mostly used 
to set up semi-classical EoS's such as the Lattimer-Swesty 
EoS [J or other modelings of the crust of neutron stars. 
We note that in the semi-classical EoS [![, the separa- 
tion between the r-cluster and the r-gas is sharp in co- 
ordinate space, but it is not so in the microscopic cal- 
culations. The idea of clusters as spatial structures in 
Wigner-Seitz cells was alrea dy p ut forth in the earliest 
studies of clusterized matter [29[ . We propose in the fol- 
lowing a method to relate this semi-classical modeling 
with the microscopic calculations where we use Woods- 
Saxon functions to fit the microscopic density in the cell. 

In the r-clustering representation, the density profile 
of the nucleons represented in Fig. [1] is decomposed in 
terms of an inside Woods-Saxon function, 



Po 



l + cxp[{r - R ws )/a}' 



surrounded by an outside gas density defined as 



ws () 

I r—qas \ J 



Pga 



1 + cxp[-(r - R ws )/c 



(14) 



(15) 



so that the sum of the two profiles gives the total profile 
p^n(r), introduced in Eq. (TTU)) . The parameters of the 
Woods-Saxon density profile are: the r-cluster bulk den- 
sity po, the radius R , the surface diffuseness a, and 
the asymptotic gas density p^f . The same quantities 
can be introduced separately for neutrons and protons, 
namely 



ws 

Pr-cl 



(r) = — 



P0,q 



exp[(r - RW s )/a^ s ] ' 



and 



p WS i 



r) 



Pgas,q 



l + exp[-(r-iS™)/a,] 



(16) 



(17) 



In Eqs. (|T6|) - (fT7|) the index q stands for neutrons (n) 
or protons (p.) The Woods-Saxon parameters entering 
Eqs. ([T4|) - ([T7| can be obtained with the help of fits on 
microscopic nucleon, neutron and proton density profiles, 
as we describe in the Appendix. An analytical model for 
the parameters will be introduced in Sec. IIVI 

Since the sum of two Woods-Saxon functions with dif- 
ferent radii and diffuseness parameters is no longer a 
Woods-Saxon function, only two of the three density pro- 
files defined above (for protons, neutrons and nucleons) 
can be employed simultaneously. In this work we con- 
sider the nucleon and proton densities as the independent 
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FIG. 1. Top: A typical microscopic density profile (triangles), corresponding to 50 protons and 1900 neutrons. The solid grey 
line corresponds to the Woods-Saxon fit of the total density, the green dashed line to the cluster density (left: r-cluster; right: 
e-cluster) and the blue short-dashed line to the gas. The density of unbound nucleons, localized or non-localized, is also shown. 
Bottom: Number of particles in volume elements as indicated, for the density profiles shown on top. The vertical lines indicate 
radial elements dr = 0.6 fm. 



functions. The nucleon and proton densities in the cell 
are defined as, 



and 



Pcell,n — P ] r-cl,n{ r ) + Pr'- a as,n{ r ) (21) 

obtained from a separate fit on the neutron density. 

On the left of Fig. [T] is shown the result of the fit for 
the nucleon density (fT8"]) in solid line, as well as the con- 
tributions of the cluster density p^- c i ( r ) an d § as density 



(18) 



pZL(r) = pV-cl, P (r) + p™_ S gas , p {r), , (19) 
respectively, while for the neutrons we choose 

pZUr) = P r {1) (r) = p ws {r) - pf s {r) (20) 
to be compared with 



pV- S gas( r ) represented respectively by the dashed and dot- 
ted lines. 

In the following we drop the label "WS" from our no- 
tations. 

The number of nucleons in the r-cluster, A r - C i, and 
the numbers of nucleons in the gas, A r - gas , are defined 
as 



A r -ci = I d 3 rp r _ c/ (r), , A r _ ffas = / d z rp r - gas (r). 



(22) 

We have similarly for the number of protons in the clus- 
ter, Z r _ c i, and in the gas, Z r 



3 5 ^r — gasi 



Z r -ci — J d 3 rp r - c i_ p (r), , Z r - gas — J d 3 rp r - gaStP (r) , 

(23) 

the latter being zero in the case of Fig. [1] Finally, the 
number of neutrons in the cluster, iV r _ c /, and in the gas, 
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N r 



are 



implying 



Nr-cl = N, 



r—cl 



r — gas 



id 



Alternatively, 



N, 



(2) 

-cl ~~ 



" ^4.r — cl Z r — C i 
Ar—gas 'Zi} — gas- 



J d 3 rp r _ cLn (r), 
Nyj gas = J d 3 rp r - gas , n {r), 



(24) 



(25) 



from a separate fit on the neutron density. 

We have checked that our numerical results, for which 
there is no proton gas, are practically the same (differ- 
ences by less than 1%) if we treat the neutron distribution 
as independent or as the difference of the total and the 
proton ones. We chose to work with the difference (i.e., 
with only two independent density distributions) to en- 
sure exact conservation of the total number of particles. 
This way we also avoid ambiguities in the definition of 
the cluster volume. 



B. Energy-space clustering 

The energy-space clustering is an alternative represen- 
tation of the density profile where clusters are defined 
as a collection of bound particles, independent of their 
localization. As a consequence, the external gas is al- 
lowed to penetrate inside the cluster. At variance with 
r-clustering, this representation allows extending the con- 
cept of clustering to supercritical densities |3(| ■ 

This representation appears natural in quantum me- 
chanics, where wave-functions with various quantum 
numbers can overlap. Since the separation between the 
cluster and the gas is described in terms of quantum 
single-particle wave-functions mainly characterized by 
their energies, as we discuss below, we refer to these two 
ensembles of particles as the e-cluster and the e-gas. 

Before we discuss the microscopic origin of e-clustering, 
we begin with an analytical parameterization of the e- 
cluster and e-gas density profiles. The density of the 
c-clustcr is defined as 



PO — Pgas 



\ + e^[{r -R ws )jc 



(26) 



and the e-gas density, Pe-^ as (r), is set to be constant in 
the coordinate variable r, 



JVS 



Pgas 



(27) 



P" J (r) 

re— gas V / 

The various parameters are the same as in the case of 
r- clustering. It is readily verified that the total density 
in the energy representation, 

pZi(r)= P Zii(r) + pl VS g a S (r), (28) 

is the same as in the case of r-clustering and yields 
Eq. (|TU1) . Therefore, the fit does not depend on the rep- 
resentation. Comparing Eqs. ((T4|) and (|26)) . we obtain 



Pgas 

Po 



i ws M 



(29) 



Ae-cl 



PO 



(30) 



where A e - C i is defined hereafter. 

Completely analogous densities p^_ d q (r) and related 
quantities can be defined for the neutron and proton 
densities separately. The neutron density profile and 
the number of neutrons in the e-cluster can be ob- 
tained in two ways, as in the case of the r-cluster, cf. 
Eqs. ([H]), (E3- Here we adopt 



N, 



- — e— cl 

N — 

Iy/ e-gas — JV e-gas ' 



A 

si-e—gas 



^e—gas ■ 



(31) 



in analogy to Eq. 

The density profiles ([26]) and (|2T[) are represented in 
the right panels of Fig. [T]by the dashed and dotted line 
respectively. The nucleon density in the cell obtained 
from the fit is represented by the solid line and is iden- 
tical to the one on the left panel. Equations (|2"9"|) . (|30p . 
(|3f [) show that there is an exact and analytical mapping 
between the two representations, and we will therefore 
use both in the following. 

To get an insight on the physical meaning of c- 
clustering, it is possible to extract the density associated 
to the bound states, the resonant states and the contin- 
uum states from the microscopic quantum calculation. 
We have classified the wave functions into three groups 
depending on their single-particle energies and space- 
extension inside the cell. First we define the quantity 
u^j as the external mean field, at a large distance from 
the cluster. A bound wave function is considered to have, 
by definition, a single-particle energy q < u^j , otherwise 
the wave-function belongs either to the resonant states 
or to the continuum states. Since in the microscopic cal- 



culation the states with e, > 



are expanded in a 



discreet-box basis, it is difficult to unambiguously dis- 
tinguish the resonant states from the continuum states. 
In our case they are differentiated according to their co- 
ordinate space extension: unbound wave-functions pop- 
ulating substantially the cell outside the cluster, with an 
amplitude comparable (within 30%) to normalized plane 
waves, are non-nocalized and are set to be in the con- 
tinuum, otherwise they are considered resonant states. 
The accuracy of this geometric method is limited by the 
wave-functions close to the centrifugal barrier which re- 
main difficult to identify. There is an ambiguity of few 
nucleons in this method which is enough for our discus- 
sion. 

On the right side of Fig. Q] we show the density pro- 
file of the continuum, or non-localized unbound states 
(squares) and that of the resonant or localized unbound 
states (dots). It is interesting to notice that the den- 
sity due to continuum states penetrates well inside the 
cluster, which confirms that the gas is not excluded from 
the central part of the cell. Additionally, one can notice 
that the continuum-state density slightly increases inside 
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the cluster compared to the asymptotic density. This 
property shall be related to the attractive nature of the 
proton-neutron interaction. The simpler approximation 
for the density of the gas inside the cluster is simply to 
set it constant and defined by the asymptotic density, as 
anticipated with Eq. ([27]) . The error introduced in the 
population of the cluster is small. 

The density of the resonant, or localized, unbound 
states is also shown in Fig. [T] There is a clear overlap 
in coordinate-space between the density of the resonant 
states and the cluster. If we want to define an c-clustcr 
in terms of its single particle states, an ambiguity exists 
concerning resonant states. These could be attributed 
to the gas, since they are unbound, or to the cluster, 
since in the e-clustering concept the gas is constituted 
by an homogeneous background, i.e., by single particle 
plane waves, while resonant states are clearly affected by 
the cluster mean-field. In order to reach a decision, we 
analyze further our microscopic calculations on a large 
number of systems in Fig. [5J The number of neutrons in 
the energy-space cluster N e _ c i is first estimated from the 
fit of the cell density and proton density via Eq. (f3~T|) . 
N e - C i is represented by the red crosses in Fig. [2] The 
number of bound and resonant states are calculated from 
the bound and resonant state densities deduced from the 
microscopic calculation as explained before. The points 
in Fig. [2] with the lower N ce u correspond to isolated nu- 
clei where the resonant states are not occupied. The 
dripping point is marked with a vertical line. Up to that 
point one observes a perfect matching between the three 
quantities represented in Fig. [2] namely N e _ c i, the num- 
ber of bound states Nb ounc i, and the sum of the bound 
and resonant states. When the first neutrons drip out 
of the clusters these three quantities become different. 
The number of states in the cluster is slightly larger than 
the number of bound states, and in general close to the 
sum of bound and resonant states, indicating that res- 
onant states should be included in the definition of the 
energy-space cluster. 

Notice that the coordinate-space representation in- 
duces naturally the idea of an excluded volume for the 
external gas, while the energy-space representation al- 
lows the cluster and the gas to overlap and there is no 
excluded volume. 

To summarize, both representations have their own 
advantages and disadvantages. They can be exactly 
mapped onto each other, and we will take advantage of 
that in the next section. However, the chosen representa- 
tion has important consequences in terms of allowed ex- 
cited states and excluded volume, meaning that all these 
aspects have to be treated consistently in the implemen- 
tation of these representations in a statistical modeling 
of supernova matter. 



IV. ANALYTICAL MODELING OF THE 
DENSITY 

In section IIII1 the coordinate-space and the energy- 
space clustering were introduced and illustrated in the 
case of 1950 Sn. In this section, we make a broader analysis 
of the results of the fit and we propose an analytical 
modeling which reproduces the parameters of the fit with 
an accuracy of the order of roughly 5% or better. 

In the analytical model the parameters of the Woods- 
Saxon function, namely the bulk density Pq S , R ws , and 
a, are expressed as functions of these variables: (A r _ c ;, 

Zr — cl: Pgasi Pgas,p) Or (^4. e _ c /, Z e ^ c i, Pgasi Pgas : p)- The 

dependence is not always direct, but goes through the 
bulk asymmetry 5 r ~ c i or <5 e _ c ;. The populations of the 
e-cluster and the r-cluster are related via Eq. (|3T)|) . Next 
we are going to explore these dependencies to complete 
the modelization of the density profiles. For this purpose 
we have performed a large number of microscopic calcu- 
lations for a selection of nuclei and cells (Ca, Ni, Zr, Sn, 
Pb) , varying the number of neutrons in the cell from sta- 
ble nuclei to approximately 3000 with a fixed size of the 
cell, with radius i? ce n=35 fm. 

A. The bulk density 

We first relate the bulk density in the clusters to the 
saturation property of nuclear matter. The saturation 
density po is defined as the density for which the binding 
energy of nuclear matter is maximized - or where sym- 
metric matter is mechanically stable, P(po,8 = 0) = 0. 
The saturation density is usually defined in symmet- 
ric nuclear matter, for which the isospin asymmetry 
5 = (N — Z)/A = 0. The saturation density can however 
be generalized to asymmetric nuclear matter as the den- 
sity for which asymmetric matter is mechanically stable: 
P(po(5),5) = 0. Apart from fluctuations arising from 
quantum shell structure, ordinary isolated nuclei aim at 
having a bulk density close to the saturation property 
of nuclear matter for an asymmetry that corresponds to 
their bulk asymmetry. In the following, it is shown that 
the bulk density in dilute nuclear clusters also tends to 
the same saturation density. 

Considering the Generalized Liquid-Drop Model 

(GLDM) [111 UK], the binding energy is expressed as a 
systematic analytical expansion around the saturation 
density po(S = 0) as, 

B(x, 6) = ^ ( c is,n + C IV , n 6 2 ) x n , (32) 

n>0 

where x = (p — po)/(3po) (here po = po{S = 0) and we 
will shorten the notation in unambiguous cases), and the 
coefficients C'is, n and Civ,n stand for the n-derivative 
of the binding energy with respect to the density p in 
the isoscalar (IS) and isovector (IV) channels. We have 
the following relations: Cis,o = B(0, 0) = Bq (binding 
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FIG. 2. For given number of protons (isotopic chain), and as a function of total amount of neutrons in the cell, shown are: 
the number of neutrons in the energy-space cluster, computed as the integral of the e-ciuster density profile and assuming a 
uniform gas, N e ~ c i; the number of neutrons in bound states, Nbound', and the number of neutrons in bound or in "resonant", 
i.e., localized unbound states, Nbound + N reaon . 



energy), Ciss = (pressure), Cisa = Koo (incomprcss- 
ibility modulus) for the isoscalar channel, and Cjy,o = Jo 
(symmetry energy), Civ,i = L (slope of the symmetry 
energy), Civ,2 = K sym (curvature of the symmetry en- 
ergy) in the isovector channel. The nuclear coefficients 
Cis,n and Civ, n give a characterization of the nuclear 
properties around saturation density. Much effort has 
therefore been invested in giving accurate experimental 
values. For instance, the binding energy is estimated 
to be B Q = -16 ± 0.5 McV, K x = 240 ± 40 MeV, 
J = 32 ± 4 McV, L = 65 ± 20 MeV [H|. For the 
SLy4 functional, in particular, the following values are re- 
ported in [H: B = -15.972 MeV, = 229.97 MeV, 
K sym = -119.74 MeV, J = 32 McV, L = 45.94 McV. 

In the following, we will write for the bulk asymmetry 
parameter 

5 r -ci = 1 — 2po.p/po 

or 

&e-cl = 1 — 2(p ,p — Pgas,p)/{P0 — Pgas) 

inside the cluster and 



in the gas. We note that the bulk asymmetry parameter 
is not equal to the total asymmetry 

Ir — cl = 1 ~ 2.Z r _ c (/y4 r _ c ; 

or 

le-cl = 1 — 2Z e _ c i/A e ^ c i 

and is more complicated to calculate, because of the pres- 
ence of a neutron or proton skin at the surface. A model 
for the bulk asymmetry parameter is presented in sec- 
tion GvH 

Setting the first density-derivative of Eq. (|32|) equal 
to zero, the saturation density in asymmetric matter 
po(5 r - c i) can be related to the nuclear coefficients Cis.n 
and Civ,n- Limiting the series expansion in Eq. (|32[) to 
n = 2, it is found indeed that 



GLDM,n=2 

Po 



(Sr-cl) = QO I 1 



ZL51 



K S 2 



, (33) 



where Qq the saturation density of symmetric matter. A 
lowest-order expression can be obtained with the limita- 
tion to n = 1 in Eq. (|32|) . 



5 r - 



gas u e—gas 



1 2 Pgas, p/ Pgas 



GLDM ,n=l I c 
Po \°r 



(34) 
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as it was already obtained in Ref. [35|, expression (64). 
Notice that, since Eqs. (pJ3"|) and (fM]) are obtained from 
considerations in nuclear matter, the density and asym- 
metry shall be those of the equivalent nuclear matter 
system, which in our notation means p r - c i and 5 r _ c ;. 

The accuracy of the expression (|33]l can be estimated 
from the comparison to microscopic calculations. The 
bulk density and asymmetry are extracted from the mi- 
croscopic calculation as the result of the fit described in 
section IIIIl The so-obtained microscopic bulk density p 
is shown in Fig.[3]as function of the respective bulk asym- 
metry S r - C i, where the various types of symbols are as- 
sociated to various isotopic chains. These are compared 
with the exact result for the SLy4 Skyrme functional and 
with the analytical expression Q33p p Q DM ' n ~ 2 (Jj r — c i), 
where the nuclear coefficients L, /Too, and K sym are set 
to be that of SLy4. The nuclear coefficients are given, 
for instance, in Table I of Ref. [HJ- The yellow band 
around the exact result expression represents a varia- 



tion of ±5%. The analytical expression p 



GLDM.n=2 



(S) 



reproduces the exact result accurately up to large asym- 
metries, while the analytical expression p^ LDM ' n ~ 1 {8) is 
only reliable at moderate asymmetries. The microscopic 
bulk densities follow the nuclear-matter saturation curve 
within a deviation of the order of 5%. Close to symmetry, 
the microscopic bulk density is systematically lower than 
Qo = po($ = 0) mostly due to the Coulomb interaction. 

In order to analyze whether the transition between 
isolated nuclei and cells plays a role, the ratio of the 
microscopic bulk density p$ over the analytical model 
Pq LDM '"~ 2 (5) are shown in Fig. 2] with respect to the 
number of nucleons inside the coordinate-space cluster 
A r - C i (top panel) and with respect to the total number 
of nucleon in the cell A to t (bottom panel). It is clear that 
the analytical model p ,n ~ 2 (S) gives an equally ac- 
curate estimation of the microscopic bulk density p for 
both isolated nuclei (A tot < 50 — 300) and nuclei embed- 
ded in the gas (A tot > 50 - 300). 



B. The bulk asymmetry parameter 

Knowledge of the neutron and proton density profiles 
requires the knowledge of the Woods-Saxon parameters 
po,q, which are related to the bulk density and the asym- 
metry parameter 6 r - c i as 



Po,n 
Po,p 



1 + 5 r -d , r s 

Po(O r -cl) 



1 — <5 r _ r 



-P0(5r-cl)- 



(35) 
(36) 



A model to estimate the isospin asymmetry 5 r ~ c i is there- 
fore needed. 

In isolated nuclei, the bulk asymmetry <5 e _ c ; = 5 r -ci 
and it is slightly different from the global asymmetry 
I e - C i = (A e _ cl - 2Z e ^ c i)/A e _ c i due to two effects: first 
the Coulomb interaction pushes protons to the surface 
of heavy nuclei, and the symmetry energy contribute to 



produce proton or neutron skins at the surface of nuclei. 
The relation between I K - r i and <5 e _ c / has been obtained 
from a liquid-drop model |36l - [38j 



Se-cl = 



Ie-cl 



3ac jZe-cl) 

SQ {Ae-ctY 



1 



9J 



(37) 



4Q (^e-ci) 1 /' 3 



where Q is the surface stiffness coefficient extracted from 
a semi-infinite nuclear matter calculation and ac the 
Coulomb parameter defined as ac = 3e 2 /(5r ). We have 
r o = 3/(47rp (0)). In Ref. |38j . a correlation pattern 
has been obtained for a large number of nuclear mod- 
els (Skyrme, Gogny, RMF) between the symmetry en- 
ergy Jo- the slope of the symmetry energy Lq and Q 
as Lq = 144.5Jo/Q — 55.5 MeV. Using this relation, the 
value for Q can be estimated for each nuclear interaction. 

The bulk asymmetry inside the r-cluster, 5 r - c i = 
1 — 2po,p/po> can be decomposed into the asymmetry 
of the e-clustcr 5 e _ c i weighted by the fraction, £ e _ c / = 
(po — Pgas)/ Po, of the e-cluster in the r-cluster plus the 
asymmetry in the e-gas 8 e ^ gas weighted by the fraction 
x e-gas = Pgas I Po of the e-gas inside the r-cluster, namely 



Pgas 



Po{$r-cl) 



Se-cl + 



Po(Sr-cl) 



(38) 



The asymmetry of the gas S gas is simply related to the 
variables of the model, 5 gas = (p ga s,n ~ Pgas, P )/ Pgas, and 
is independent of the representation. The bulk density is 
given by expression (|33j) . and the cluster asymmetry 5 e -ci 
is obtained from Eq. ([3"T]) . Notice that po is a function of 
S r - C i and therefore Eq. (|38|) is a self-consistent equation, 
to be solved by iteration as we will discuss later on. 

The bulk proton fraction y p , r -ci is obtained from the 
asymmetry parameter 5 r - c i as y p . r -ci = (1 — 5 r - c i)/2. 
The ratio of the proton fraction y Pt r-ci obtained from 
the analytical model over the microscopic one is shown in 
Fig. [5] with respect to the total number of neutrons N to t 
shifted by the number of neutrons at the drip-line for 
each isotopic chain NdHp (top panel), and with respect 
to the number of neutrons inside the coordinate-space 
cluster N r - C i shifted by Nd r i P - The results justify our 
approach. A deviation observed arount the drip point is 
cured as more neutrons are added to the gas. 

One might consider to apply expression (|37| to the r- 
clustcr instead. In such a case, we found that, compared 
to the microscopic calculations, the asymmetry is over- 
estimated as the density of the neutron gas increases. 



C. The surface diffuseness 

In ordinary nuclei, the surface diffuseness a is of the 
order of 0.55 ± 0.05 fm [3^| and this value is different for 
light nuclei and heavier stable nuclei. There is therefore 
a dependence on the mass of the nuclei. We now turn 
to show that this mass dependence can be largely ab- 
sorbed in a dependence on isospin, both in the case of 
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FIG. 3. Central (bulk) density of clusters from the fits of microscopic density profiles (points), compared with the saturation 
density of nuclear matter at a given asymetry, given by: an exact numerical calculation, solid red line; approximation (|33l) . 
short-dashed blue line; and approximation (j34[) . dotted gray line. 
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FIG. 4. Comparison between the values of the bulk density as obtained from the fits and from the saturation property of 
nuclear matter, as a function of the r-cluster population or the total amount of nucleons in the cell. 



isolated nuclei and nuclei embedded in a gas [40J. The 
set of nuclear systems studied in this work includes a 
few neutron-deficient nuclei, stable nuclei and neutron- 
rich nuclei or in-medium clusters. In Fig [6] we show the 
surface diffuseness obtained from the fit to microscopic 
density profiles for the same set of nuclear systems as in 
Fig. [3] The diffuseness of the nucleon density is shown 
on the left panel. 

An almost linear correlation between the surface dif- 
fuseness a n and the bulk asymmetry (5 r - c i) 2 is apparent 
and independent from whether we consider isolated nu- 
clei (colored shapes) or clusters embedded in a gas (black 
shapes in gray disks). It has therefore been fitted as, 

a = a + /36t_ cl [fin]. (39) 

The parameters a and j3 are given in Table [I] consider- 
ing various sets of nuclei in the fit: Ca, Ni, Zr, Sn, Pb, 
and finally taking into account all the nuclei in the set. 



Linear fits were performed for the proton and neutron 
diffuseness parameters too, 

a p = a p + (i p 5l_ cl ; a n = a n + /3 n 6^_ d [fm] (40) 

and the results are also given in Table [I] In all cases the 
standard deviation 




N 

5>-(a + /3,5 r 2 _ ci )P, 



where N the number of points in each sample and 
the diffuseness, at each point, from the microscopic fit, 
is equal to about 0.2 for a p and 0.3 for a, a,,. If the fit is 
performed by excluding N < Z nuclei from the sample, 
very similar results are obtained, namely a = 0.53, (3 = 
1.06, a p = 0.53, P P = 0.35, a n = 0.54, /3„ = 1.14. 

The various parameters are very stable and depend 
only weakly on the chosen set of input. This indicates 
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FIG. 5. Quality of our approximation for the proton fraction in the interior of the cluster - equivalently, the bulk asymmetry. 
The ratio of the proton fraction, computed via the analytical model, and the value fitted to the microscopic calculations, as a 
function of: the neutron excess or deficiency in the cell with respect to the drip line (left); and the neutron excess or deficiency 
in the r-cluster (right). 




FIG. 6. Linear fit of the diffuseness parameter corresponding to: the nucleon density profile (left); the proton density profile 
(right-top); and the neutron density profile (right-bottom). Colored shapes correspond to isolated clusters, while black shapes 
in gray disks correspond to clusters in the presence of a gas. 



that a, /3 and the proton and neutron counterparts might 
be related to general properties of the nuclear interaction. 
For comparison, we have performed a fit also for a simi- 
lar set of calculations using the LNS interaction [4l| • We 
found in that case a = 0.53 and /3 = 0.98, i.e., a some- 
what smaller diffuseness, especially for the neutron-rich 
systems. However, these values are still compatible with 
the ones obtained with SLy4, considering the scattering 
of pseudo-data points, quantified by a. All values are 
compatible also with the values extracted from experi- 
mental data. 

We notice from the values discussed above and from 
the right panels of Fig. [5] that the diffusness parameters 
of the proton and neutron density profiles arc different, in 
concordance with a recent study of the difference between 



the proton and neutron surface thickness [38j . The two 
possible reasons are the presence of a neutron or proton 
skin at large asymmetries and the Coulomb interaction 
spoiling the symmetry between the species. The proton 
diffuseness a p is in general lower than the neutron diffuse- 
ness a n , except for the smaller bulk asymmetries, which 
correspond to neutron-deficient nuclei. The difference in- 
creases for neutron-rich systems. We have checked that 
the strongest effect is due to the nuclear interaction. We 
have performed calculations and fits by switching off the 
Coulomb interaction and obtained very similar results, as 
shown in Fig. [6l This means that in proton-rich systems 
(negative bulk asymmetry) the role of the protons and 
the neutrons would be roughly reversed, i.e., the proton 
diffuseness parameter would be larger. 
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TABLE I. Parameters in Eq. (|39l) from linear fits to various sets of nuclear systems: Ca, Ni, Zr, Sn, Pb and from all sets. 



1.1 

1.05 
1 

0.95 
0.9 

0.85 
0.8 





'l-(jt 2 /3)x' — 




N=Z * 




Ca * 




Ni 




Zr □ 




Sn ■ 

Pb ° 







0.005 0.01 0.015 0.02 0.025 0.03 

(a ffi /R HS,fi t) 2 



FIG. 7. Quality of the relation flJTJ. 



Note that for symmetric systems, and ignoring the ef- 
fect of the Coulomb interaction, the diffuseness param- 
eter should not depend on the nucleon species. All our 
results are consistent with a rounded value a = a p = 
a n = 0.53. A dedicated investigation of the dependancc 
of all the parameters on the interaction properties would 
be worth undertaking in the future. 



D. The radius of the cluster 

The radius parameter R ws entering the density profile 
is defined as: 



R 



\vs 



R 



HS 



1 



7T 2 / a N 2 



3 \R HS 



(41) 



with the help of the equivalent homogeneous-sphere 
value, which is given by: 



R 



HS 



HS 



\ 47T 



1/3 



(42) 



The equivalent homogeneous-sphere volume, V HS , is 
simply expressed as 



V 



HS 



A r —ci 
Po(S r -c 



= V d 



(43) 



and defines the volume of the cluster, V c i. 

Eq. (|41j) is a series expansion in the small parame- 
ter a/R HS . In order to check the accuracy of Eq. (f4"Tj) 
the ratio R ws /R HS is shown in Fig. [7] as a function of 



(a/R HS ) 2 and for the whole set of nuclei considered in 
this work. The various symbols correspond to the iso- 
topic chains for Ca, Ni, Zr, Sn, and Pb, while the expres- 
sion f4"Tj) is represented by the solid line. The matching 
between the symbols and the solid line justify the trun- 
cation in Eq. (|4T|) . 

A parameterization of the radius of the Woods-Saxon 
profile is obtained by injecting expressions (f3"3")l . (|42p . 
and (|4"3"|) inside Eq. (|4"Tj) . This model is compared to 
the Woods-Saxon radius in Fig. [8] where the ratio of the 
analytical radius (j4"T|) over the Woods-Saxon radius from 
the microscopic calculations is represented. For the com- 
pleteness of the analysis, the ratio is represented with re- 
spect to various quantities: 5 r - c i, A r - C i, and A to t- The 
same set of nuclei as in Fig. [3] is taken into account. The 
yellow band give a deviation of ±5% around the analyt- 
ical expression. The analytical radius f4"2")l gives a good 
estimation of the Woods-Saxon radius. 



E. The analytical model in practice 

We can now summarize the basic equations in the 
model and how they are used in practice. The nu- 
cleon density profile is determined by Eqs. ([T3 ]) .([I5 ]) .(fT5 ]) 
or, equivalently, Eqs. ([2rj |) .([27 ]> .(|2g ]) . The variables on 
which the Woods-Saxon parameters depend are chosen 
to be the composition of the cluster, which depends 
on the representation, and the gas density, which does 
not. In particular, the variables are either (A r - C i, Z r — C i, 
Pgas, Pgas, P ) for the r-cluster, or (A e _ ci , Z e _ d , p gas , 
Pgas.p) for the e-cluster. Equivalently, the gas asymmetry 
Sgas = 1 - 2pgas,p/Pgas may replace Pgas, P as a variable. 
The corresponding density profiles are given by Eqs. (|15[) 
and (|26[) . respectively. For the various parameters we 
now have: 



po(5 r - c i) = go 1 - 



3L 



K t) 2 



a = a + /35 r _ c i, 



with 



R ws = R HS 



R HS = 



7r 2 ( a 



3 \R HS 



3A 



r—cl 



47rp (£ 



r—cl j 



1/3 



(44) 
(45) 
(46) 

(47) 
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FIG. 8. Comparison between the cluster radius determined through the analytical model and the radius obtained from the 
microscopic density profiles, as a function of relevant quantities, namely: the bulk asymmetry of the r-cluster; the population 
of the r-cluster; and the total number of particles in the cell. 



The r-cluster asymmetry entering the above expressions 
is obtained self-consistently with the help of Eq. arL d 
the relations 



K.-rl = 1 



Po($r-cl) 

Pgas,p 
p0,p(Sr-cl) 



Ar—cl'i 
Z r —cl'i 



Ie-cl + 



3ac (Z B - C ,y 



'(A e _ e ,)V3 

1 4. 9J o 1 

1 + 4Q (A c _ ci ) 1/3 



(48) 

(49) 
(50) 



Po{$r- 



S r . cl = ( 1 - -#!£_ ) Se- d + -£^<W (51) 



Po(S r -d) 



The procedure depends on the chosen variables. If 

— cl i ^r — cl j Pgas i Pgas,p 

) are given, we proceed as 
follows. We start with the initial estimate <5, c i(o) = 

1 — 2Z r _ c ;/A r _ c /; we use Eq. (|44| to estimate p^; then 
Eqs. (gSJ), ((33) estimate the e-cluster population and next 
Eq. (|50p its asymmetry; Eq. (f5"Tj) gives a new estimate for 
the asymmetry, 5 r _ c ;n); and so on until convergence. If 

i Pgas i Pgas,p 

) are given, one readily obtains 
8 e -ci from Eq. (|50p and can proceed as follows. An initial 
estimate for pjf is the saturation density of symmetric 
matter, p^ = qq. The r-cluster asymmetry is then esti- 
mated from Eq. (|5"Tj) . This leads to a new estimate Pq 
via Eq. P4")l . And so on until convergence. 

The parameters determining the proton and neutron 
density profiles are readily obtained. Proton and neu- 
tron bulk densities are determined by the total bulk den- 
sity and the bulk asymmetry Along with the numbers 
of neutrons and protons in the cluster they yield the re- 
spective radii. Diffuseness parameters are given by the 
approximation P0|) . with a PiTi and /3 PtTl possibly depend- 
ing on the interaction, but in general close to the values 
given in Table HI as already discussed. Note that for 
proton-rich systems (S r - C i < 0) the proton and neutron 
diffuseness parameters are exchanged with respect to the 



values given in the Table. 

In applications relying on an excluded volume, care 
must be taken to define the volume of the cluster (both 
protons and neutrons) uniquely. A straightforward strat- 
egy is to set the proton and neutron r-gas inner radius 
equal to the total nucleon r-gas inner radius. The proton 
and neutron density profiles would still be determined by 
the above relations, but a reshuffling of nuclcons between 
the gas and the cluster population might be necessary. In 
our numerical calculations this ambiguity is avoided be- 
cause there is no proton gas and we do not treat the neu- 
tron density profile as independent, see Eqs. ([20"]) . (f23|. 



V. ENERGY OF A CLUSTER IN A GAS 

We now turn to the energies of clusters embedded in a 
nucleon gas. We will define these energies based on our 
two representations for dilute clusters and then study 
them with the help of the analytic density profiles and 
the local density approximation (LDA). 

The LDA allows us to go beyond the restricted do- 
main of Hartree-Fock, namely to study combinations of 
cluster and gas compositions that would not be gen- 
erated through energy minimization. It is well known 
from present EoS models that a very large distribution 
of clusters of any size and isospin is expected in the dif- 
ferent thermodynamic conditions relevant for supernova 
physics. Information on what are essentially excited con- 
figurations is thus mandatory for the treatment of clus- 
terized matter at finite temperature. In fact, we will 
show that relying solely on Hartree-Fock or, in general, 
variational calculations (similarly: existing nuclei) can be 
highly misleading when attempting to generalize the re- 
sults to more exotic systems. Such results remain of value 
at certain domains of low temperature, but, in general, 
one should keep in mind the associated limitations when 
extrapolating. 

With practical applications in mind, in particular NSE 
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implementations [10( , where the energetics of the gas and 
of the cluster (whether in energy or coordinate repre- 
sentation) arc treated separately, we will partition the 
total nuclear energy of the Wigner-Seitz cell into the 
nuclear energy of the cluster and that of the gas. In 
particular, the energy due to the interaction of the two 
coexisting systems will be conveniently assigned to the 
cluster. The gas can then be treated as a uniform sys- 
tem - for example, as is currently done in the statis- 
tical model of Refs. [HI, within the temperature- 
dependent Hartree-Fock method. Then we need only be 
concerned with the energetics of the cluster. 

Our first task is therefore to disentangle the energy 
of the cluster from that of the gas. For that, let 
us call E ce u the nuclear energy of the whole system 
(A r _ cl , Z r _ ch p gas ,6 gas , R ceU ), as calculated from a mi- 
croscopic approach. Our aim being to describe arbitrary 
configurations, that is excited states as well as ground 
states, this is given by the LDA approach in terms of 
the analytical Woods-Saxon density profiles, developed 
above. The total energy can be partitioned, without in- 
troducing any ambiguity, as the sum of the bulk energy 
of the r-gas, the bulk energy of the r-cluster, and a cor- 
rection SE, which we will identify as the surface energy: 

Ecell = (V C ell — V c l)e(pgas,8 gas ) 

+V c ie(p ,5 r - Cl ) + 5E. (52) 

Here e{p, S) is the energy density of homogeneous matter 
at density p and asymmetry d. In the Skyrme model it 
equals the energy density of Eq. ©, with the gradient 
terms vanishing. In the above, V c i = ^ttR 3 hs is the vol- 
ume of the cluster and V ce ii = ^R-ceii ^ s * nc v °l urne of 
the cell. 

As already argued, we may absorb SE into the energy 
of the cluster, E r _ c i, and henceforth be concerned only 
with this quantity, i.e., 

E r -cl = E ce u — (V C ell — V c l)e(pgas,Sgas)- (53) 

For Skyrme functionals within the LDA, the above yields 
readily an integro-analytical expression for the energy of 
a cluster (A, Z) embedded in a gas of given density and 
asymmetry: 

E<r—cl(A r — c l : Z r — C i 7 Pg as , 5g as *) 

4tt [ R ~ U e(pY e i S i(r),5 ce ii(r))r 2 dr 
Jo 

-H R lell 

-[R HS {A r - c i,Z r - ch p gas , 5 gas )] 3 }e(p gas , S gas ), (54) 

where the local density p^f as well as the local asym- 
metry 

5 C eii(r) = l-2p% l l p (r)/pZ$(r) (55) 

(see also Eqs. (fT5|) . (fT5j) ) are determined by 
(A r - C i, Z r - C i, p gas ,6 gas ) within our analytical model. 



Quite similarly, we readily have an expression for the 
energy of the e-cluster: 

Ee — cl(,A e — c i 7 Z e — c i 7 p gaS7 5 gas *) 

47r / P^eU ( r 5 A e - c i,Z e - c i,p gaS7 5 gas )r 2 dr 
Jo 

-^Rl ell e{p gasi 5 gas ), (56) 

where the cluster-gas interaction energy in the region 
where they coexist, namely V c i, has been absorbed into 
the energetics of the cluster. There is no ambiguity in 
this choice, since the total energy in the cell is conserved 
by construction. In the above expressions, the radius of 
the cell Rceii , chosen randomly and much larger than the 
cluster's radius, is relevant only to the gas energy and 
therefore is an auxiliary quantity, not affecting the ener- 
getics of the cluster, in either representation. 

Next we will discuss our numerical results based on 
the above results. First, we will discuss the total cluster 
energy in order to validate the LDA. Next, we will discuss 
separately the bulk and surface energy of the r-clustcr. 

A. Total cluster energy 

In Fig. [9] wc show the nuclear part of the cluster en- 
ergies resulting from the Hartree-Fock calculation, as 
points connected with lines. They are computed using 
Eq. (|53[) with E ce a the Hartree-Fock energy of the cell. 
We notice that the binding energies of the clusters de- 
crease as the number of nuclcons in the system increases. 
Indeed, for a given isotopic chain, A to t controls the to- 
tal number of neutrons. The larger its value, the more 
asymmetric is the cluster that is favored by energy mini- 
mization, which is reflected in the decreased binding. We 
will return to the role of the asymmetry when we discuss 
separately the volume and surface energy. 

We now ask whether the LDA is an acceptable approx- 
imation to the microscopic energy. Therefore, we will 
compare the Hartree-Fock energies of the clusters with 
their LDA energies, computed from Eq. (|53")) with E ce u 
being replaced by the LDA energy of the cell, using the 
analytical density profiles - equivalently, Eq. (fMf. The 
difference of the two energies per particles is also plotted 
in Fig. |H1 It corresponds to the points just above zero in 
the figure. We notice that the difference remains roughly 
the same for all isotopic chains and A to t (or asymme- 
tries), and regardless of the presence or not of a gas, at 
the level of half an MeV per cluster particle. This is a 
very useful outcome: if we are to quantify additive in- 
medium modifications to nuclear energies, the LDA will 
be as reliable as Hartree-Fock at the very least. 

B. Bulk energy 

Even within the analytical density profiles in the LDA 
approximation, the in-medium modified cluster energy 
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FIG. 9. Points connected with lines: Energy of r-clusters, 
determined by subtracting the r-gas contribution from the 
total Hartree-Fock energy in the cell. The points in the upper 
part indicate the deviations when the total energy of the cell 
is computed within the LDA and using the model density 
profiles. 



functional is a complex integral function of four variables, 
namely the cluster atomic and mass number Z and A, as 
well as the density p gas and asymmetry S gas of the sur- 
rounding gas. To simplify this expression and get an 
insight into its physical meaning, we will distinguish be- 
tween the volume energy and the surface energy of the 
cluster. A separation into surface and volume terms will 
also help develop formal extentions of the liquid drop 
model, in the form of, e.g., Sai(p gas ,S gas ), where at a 
liquid-drop parameter and Sat its correction, depending 
on the composition of the gas at the very least. Such ad- 
ditive corrections could be plugged into any mass table 
or model in an intuitively simple manner. 

The bulk (or volume) energy per particle of the r- 
cluster is defined as 



■E<vol\Ar—eli &r—cli Pgasi &gaB) — 



-e(po(S r - c i), S r - C i). 



Po(8r-cl) 

(57) 

The bulk asymmetry of the cluster S r - C i, correspond- 
ing to given variables (A r _ c i, Z r _ cl , p gas ,8 gas ), can be 
computed through our analytical model and determines 
uniquely the bulk density as well. The calculation of the 
bulk energy with the help of the Skyrme energy density 
is then straightforward. The result is a function of the 
bulk asymmetry, and is shown as the solid line in Fig. 1101 
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FIG. 10. The lines show the volume energy per particle in the 
r-cluster, for the SLy4 functional, according to our model (red 
solid line) and according to a liquid-drop model, quadratic in 
Sr-ci (dotted line). The points show, for the various isotopic 
chains studied in this work, the deviation of our model volume 
energy per particle from the values obtained directly from the 
Hartree-Fock calculations. 



A quadratic function, corresponding to a simple liquid- 
drop model for the volume energy and with parameters 
which are consistent with the SLy4 functional [43| , is also 
plotted (dotted line). At low asymmetries the two curves 
coincide. The role of higher-order terms in the asymme- 
try, including those due to the kinetic energy, becomes 
apparent at large asymmetries. 

We take the opportunity to pcrfom one more test of our 
analytical modeling of bulk cluster densities. In Fig. [TU] 
the points correspond to the cluster configurations which 
we obtained within Hartree-Fock. They show the devia- 
tion between a) the bulk energy per particle, calculated 
through the analytical density model as described above, 
and b) the microscopic bulk energy per particle, namely 
the one calculated using those values for the bulk density 
and asymmetry, which resulted from a fit to the micro- 
scopic density profiles. We notice that the deviations 
are consistent with zero and in general there are no sys- 
tematic deviations as we go to large asymmetries. Small 
deviations corresponding roughly to the drip point are 
cured at higher asymmetries (cf Sec. II V Bl and Fig. [5]). 



C. Surface energy 

We now proceed to a study of the quantity 5E, which 
is the difference of the total energy and the bulk energy 
and which wc will identify as the surface energy. We will 
study this quantity in terms of the r-cluster. Switching 
to e-clustcr variables is straightforward, thanks to the 
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geometric relations between the two types of clusters. 

In principle, SE depends in unknown ways on the vari- 
ables (A r - cl , Z r - C i, p gas ,5 gas ). If it is a proper surface 

quantity, however, it should scale as A 2 J^ cl . We further- 
more expect that <5 r _ c ; can replace Z r _ c \ as the relevant 
cluster parameter. We are then left with the variables 
(S r - c i,p gas ,S gas ). Regarding the quantity 



s[5 r — cli Pgasi ^gas) 



SE 



A 



2/3 



the following limiting case should hold: 

s(S r -cl, Po(Sr-cl), Sr-cl) = 0, 



(58) 



(59) 



because in the case where the cluster and the gas have the 
same density and asymmetry we arc in the limit of homo- 
geneous matter and therefore the surface energy should 
vanish. Furthermore, for relatively small absolute asym- 
metries, 

s(S r -cl,Pgas,Sgas) ~ s{—6 r - d , p gaa) —6 gas ). (60) 

The weak violation of this relation is due to Coulomb 
effects. 

We will now focus on a neutron gas, S gas = 1, though 
very similar studies can be performed in other cases. In 
short, we are interested in the quantity 



S[0, — cl-, Pgas) — ^2/3 



(61) 



r — cl 



assuming a neutron gas. 

In principle, the functional dependence of s on the two 
variables is not unique. Different clusters (A r _ c i, Z r ^ c i) 
may correspond to the same combination of (<5 r _ c ;, p gas )- 
The uniqueness will have to be demonstrated numeri- 
cally. To this end we proceed as follows: For given p gas 
we obtain the analytical density profiles corresponding to 
nuclei with many different combinations of (A r - C i, Z r - C i) 
and with a large variation in particle number, namely 
from 50 to 400 particles. This provides us with a value 
for the bulk asymmetry. We also calculate SE as al- 
ready demonstrated. We thus obtain, for given 5 r ^ c i and 
p g as, a number of surface-energy results corresponding 
to clusters of different populations. These "data points" 
are displayed in Figs. [TT1 and [121 and discussed below. 

Let us begin with Fig. [TTJ where s(S r ^ c i 1 p gas ) is plot- 
ted as a function of 5 r - c i for the indicated values of gas 
density (given in fm~ 3 ). For low gas densities the points 
form well defined lines, and even at higher gas densi- 
ties the dispersion of points is moderate. We conclude 
that the scaling with respect to A^,_ , indeed holds to 
a satisfactory degree. This important result means that 
the definition we have introduced, eq.([5*T|). correctly rep- 
resents the bulk part of the energy of the cluster. We 
notice that eq. ([57)) has no dependence on the external 
gas. This means that the presence of a medium does 
not influence the bulk term, which keeps being defined 
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FIG. 11. Model surface energy for the SLy4 functional, ob- 
tained as described in the text, for given neutron-gas den- 
sity (values indicated in fm -3 ) as a function of bulk r-cluster 
asymmetry. The secondary line for p gaa — 0.02 fm~ 3 , marked 
(0.02), corresponds to the ansatz ( 1031) . 



by the local nuclear matter saturation condition what- 
ever the external medium. This demonstrates that the 
energy increase of the cluster with increasing neutron gas 
observed in Fig. [5] is not an in- medium effect, but simply 
reflects the decreasing saturation density (and increas- 
ing surface energy, see below) associated to an increasing 
bulk asymmetry. 

The points forming the line for p gas = 0.0 correspond 
to isolated nuclei. For small asymmetries, a roughly 
quadratic law is observed. Interestingly, the line is con- 
vex around zero asymmetry, implying, in a liquid-drop 
picture, a positive surface-symmetry coefficient, in quali- 
tative a gre ement with certain calculations in semi-infinite 
matter [44j, but in contrast with many evaluations from 
the literature 4^, |46| , including ones based on the SLy4 
interaction [34 1. Our result, which we found quite robust 
with respect to approximations involved (for example, 
the use of LDA instead of Hartrcc-Fock energies), war- 
rants further investigations in the future. To some extent 
the discrepancies may be due to employing different def- 
initions for the nuclear surface and different partitions 
between the volume and surface energies [13, [4(| , as 
well as due to our use of the bulk asymmetry S rather 
than the global asymmetry / as the relevant parameter. 

We observe that, as the neutron-gas density increases, 
the surface energy decreases. On the neutron-rich side 
this is easy to understand if we consider the limiting case 
of a gas with the same density and asymmetry as the 
cluster, defined by Eq. (|59[) . Furthermore, we observe 
that the various lines are not symmetric around zero. 
This should be expected, since a proton-rich cluster in 
a neutron-rich gas is physically a very different entity 
from a neutron-rich cluster in a neutron-rich gas. It 
is worth noting that the functional dependence at low 
or negative asymmetries could not have been inferred 
from our Hartrce-Fock calculations, which naturally favor 
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FIG. 12. Thick labelled bands of points: Model surface energy 
for the SLy4 functional, obtained as described in the text, 
for given bulk r-cluster asymmetry (values indicated) as a 
function of neutron-gas density. Thinner lines: Ansatz (|64l) . 
for bulk asymetries equal to 0, 0.2 and 0.4. 



neutron-rich clusters. Only the analytical model for the 
density profiles could probe this extended domain. The 
only assumption related to energy minimization is that 
the bulk density of the cluster equals the saturation den- 
sity at the given value of asymmetry. A different bulk en- 
ergy would imply compression (or decompression) , which 
would cost too much energy to be counter-balanced. 

Of course, symmetric (not to mention, proton-rich) 
clusters would not survive in a neutron gas in most sit- 
uations of interest, in particular low-temperature asym- 
metric matter in the crust of a neutron star. At suffi- 
ciently high temperatures, though, such configurations 
could become active and their correct treatment might 
be important. 

We now turn to Fig. 1121 where the quantity 
s(S r - c i, p g as) = &E/A 2 J^ cl is plotted for the indicated val- 
ues of asymmetry as a function of gas density. We notice 

2/3 

again that it scales well with A r _ cl and it generally de- 
creases for higher gas densities. We stress once more that 
the functional forms for the lower asymmetries would not 
have been possible to obtain based on an unconstrained 
variational calculation. 



The functional form of the surface energy obtained 
here goes beyond the somewhat heuristic approaches to 
the surface energy introduced in Rcfs. 

0,113 • As a gen- 
eral statement, it appears that simple analytic estima- 
tions of surface corrections can be largely inaccurate and 
the complete result from the LDA has to be considered. 
As an example, one may suppose that the in-medium 
corrections can be simply taken into account through the 
transformation from r-clusters to e-clusters. Indeed, we 
have seen that e-clusters can be interpreted as the en- 
semble of single-particle wave functions (bound and res- 
onant) in the cluster region, with the exclusion of the 
continuum states that represent the external gas. In this 



picture one should have: 

5E = a s (5 e _ cl )A 2 e L 3 cl , 



(62) 



where a s can depend on the cluster asymmetry but 
should not depend on the gas. Since A e _ c ; = (1 — 
Pgas/ Po(Sr-ci))A r - c i, the uniform- matter limit (JSHJ) is 
automatically satisfied. Replacing into eq. (|62)) we get 



S\$r—cli Pgas) — 



1 - 



Pgas 



Po(5 r -c 



2/3 



S(d r - cl ,0). (63) 



The results for s(8 r — c i, Pgas = 0.02fm~ 3 ) are shown in 
Fig. HU (line marked "(0.02)"). For extremely neutron- 
rich clusters the ansatz works quite well, but, being sym- 
metric around zero asymmetry, it fails to reproduce the 
correct trends at lower asymmetries. 

Finally, it is worth commenting on an ansatz which we 
tried at earlier stages of this work J48|. Inspired by the 
liquid-drop formula and the conditions (|59p . (|60p . it was 
proposed: 



s(S r -cl, Pgas) = ( 1 



Po{8 r - c i) 



[a s + a ssym 8 2 _ cl \ 



a. 1 - 



Pgas 



Po(Sr-cl) 



~\~&ssym ( &r — cl 



Pgas 



Poi&r-cl) 



(64) 



where a, and 



l ssym 



liquid-drop constants such that the 



surface energies for p gas = are reproduced. The 
surface-symmetry parameter a ssym was therefore posi- 
tive. An excellent agreement of this ansatz with the 
Hartree-Fock results of neutron-rich clusters in a neutron 
gas was observed [48| . We now notice from its analytical 
form, that this ansatz predicts a larger surface energy 
for a proton-rich cluster than for a neutron-rich one (in a 
neutron gas), i.e., the opposite effect than what our LDA 
approach gave. This is not to say that proton-rich clus- 
ters in a neutron gas are a particularly relevant config- 
uration, however this outcome demonstrates the risks of 
relying on Hartree-Fock results to derive general trends. 
The disagreement between our ansatz and the LDA re- 
sults is further demonstrated in Fig. [12] The surface 
energy given by Eq. (|64p with a s = 21.5 and a ssym = 32 
corresponds to the thinner curves (the scaling with A 2 / 3 
is respected by construction) in Fig. [T^J The choice of a s 
and a ssym is such, that for isolated nuclei, p gas = 0, the 
LDA and the ansatz give the same result. For low gas 
densities, below 0.02fm~ 3 , there is reasonable agreement. 
It is obvious, however, that there are large deviations in 
various domains of gas density and asymmetry. 



VI. CONCLUSIONS 

This work presents a quasi-analytic modelling of the 
in-medium modifications to cluster density profiles and 
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energies, when these latter are immersed in a dilute nu- 
cleon gas, as it is the case for the stellar matter which 
is produced in the core of supernova and in the crust of 
neutron stars. 

Our approach is valid in principle for all but the light- 
est nuclei, where there are not enough nucleons to yield 
saturated matter. The density profile in a wide num- 
ber of nuclei and Wigncr-Scitz cells, as obtained from 
Hartree-Fock calculations with the SLy4 effective inter- 
action, have been analyzed and fitted through Woods- 
Saxon functions showing a continuous description from 
isolated nuclei to nuclei in a medium. 

We have discussed two alternative possible representa- 
tions of clusters inside a medium, depending on the local 
density (coordinate-space clusters) or on the energy of 
the single-particle states (energy-space clusters). Using 
these two representation and the mapping between them, 
we have been able to propose an analytical self-consistent 
model to relate the parameters of the Woods-Saxon func- 
tions to the global variables of the EoS. 

The energies of clusters in a neutron-rich gas were 
studied with the help of the local-density approxima- 
tion, which allowed us to investigate configurations which 
are not explored by the Hartrcc-Fock (in general, varia- 
tional) calculations, but are of importance in a finite- 
temperature stellar environment. 

The quality of the LDA approximation was tested on 
the Hartree-Fock sample. Not surprisingly, this approxi- 
mation systematically overestimated the microscopic re- 
sult. However, the in-medium cluster energy shift is re- 
markably well reproduced, showing that the influence of 
the medium is essentially due to the modification of the 
local density profile. We have discussed in detail the in- 
medium effects due to a pure neutron gas, which is rel- 
evant for the physics of the neutron star crust. For this 
case, we have demonstrated through the comparison with 
the microscopic calculations that neither the isoscalar nor 
the isovector part of the bulk energy of the clusters in the 
coordinate space representation are affected by the exter- 
nal gas, and the global bulk energy only depends on the 
bulk asymmetry of the cluster. As a consequence, the 
excluded-volume approach appears as a reasonable zero 
order approximation to account for bulk in-medium ef- 
fects. An important modification of the surface energy 
of the clusters is however observed. A decreasing surface 
tension with increasing density of the gas is observed for 
all cluster asymmetries. This leads to a sizeable modifica- 
tion of the cluster energy functional in a dense medium, 
which will have to be accounted for to have reliable EoS's 
for astrophysical applications in a near future. 

In the future, we shall add the Coulomb corrections 
to the proton mean-field as well as the contribution of 
the pairing correlations to the density profiles and the 
Wigncr-Scitz binding energies. 

Acknowledgements: Work supported by the ANR 
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Appendix A: Woods-Saxon fit of the microscopic 
density 

Here we describe how we fit a microscopic density 
profile, corresponding to nucleons, protons, or neutrons, 
with a Woods-Saxon analytical profile, 



PO ~ Pgas 



l + cxp[(r-i?^ s )/a] 

PO 



Pgas 



Pga 



1 + exp [(r - R ws )/a] 1 + cxp [-(r - R ws yffl 



(Al) 



The microscopic density p™/^^) is given on a radial 
mesh of N points with = iR s t ep spanning a Wigner- 
Seitz cell of radius R ce ii = NR s t ep - During the fit we 
wish to conserve the total number of particles, 
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(A3) 

within a given tolerence 5 A of, typically, half a particle, 
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(A4) 

The optimal profile, within this restriction, shall be de- 
termined by minimizing the standard deviation 
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The parameters to be determined arc the bulk and gas 
densities po and p gas , the radius R ws and the diffuseness 
parameter a. At the beginning of the fitting procedure 
the standard deviation is initialized at a large number. 

The radius corresponds to the inflection point of the 
Woods-Saxon profile, 



±r0 



WS, \ 
U V) 







(A6) 



and can be readily obtained by numerical differentiation 
of the microscopic density. All radial derivatives are cal- 
culated at 5-point precision. The radius will lie between 
two points for which the second derivative changes sign. 
Its precise value is obtained by interpolating between the 
two points. We then determine the value of the density 
as well as its derivative at this point, by interpolation, 
P™u( rWS ) and d P"eu( r )/ dr \r=Rws, respectively The 
value of the radius as well as the density and its deriva- 
tives are thus determined directly from the microscopic 
density profile and remain fixed during the fit. 

Next we shall vary po an d p ga s within limiting val- 
ues, [po,min,PO,max] and [p gas ,min, Pgas,ma,x] respectively, 

determined from the microscopic profile. The micro- 
scopic density outside the cluster is not perfectly homoge- 
neous, but shows weak oscillations and more so towards 
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the radius of the Wigner-Seitz cell. The oscillations are 
largely an artefact of the microscopic calculation and the 
imposed boundary conditions. We therefore determine 
the limiting acceptable values for p gas as the minimum 
and maximum numerical value of p mlc (r) within a radius 
[2R ws ,R cell ], Pga s,min and p ff as, max, respectively. 

A typical density profile, determined through Hartrcc- 
Fock, will show ripples in the interior of a nucleus or 
cluster. Initially we determine the minimum and maxi- 
mum numerical value of p™n i r ) within a radius equal to 
0.7R WS , Po,min,o.7fi ws and Po,max,o.7fl ws > respectively. If 
we determine the minimal value of the density up to a 
value too close to R ws , we run the risk of accidentally 
accepting a too- low value of po. On the other hand, we 
wish to avoid that the minimal and maximal values are 
too-restricted, due to their limitation up to 0.7R WS . In- 
deed, it is likely that the shoulder of the density profile 
is critical for a correct determination of po- Therefore 
we extend the allowed interval by considering that, for 
genuine Woods-Saxon profiles, we have: 

pZl S l (R WS ) = \(P0 + P 9 as) => PO = 2PZUR WS ) - P g as- 

(A7) 

Taking this into account, the minumum and maximum 



acceptable values for po are determined by 

P0,min = min{p min o , 7R WS, 2p mlc (R W ) - p ga s, max} 
PO.max = niax{po,max,0.7i? ws , 2p mlc (R WS ) — p g as,min} 

(A8) 

Finally we vary, within the acceptable intervals and by 
sufficiently small steps, both po an d p ga s in nested loops. 
At each step we proceed as follows: the diffuscness pa- 
rameter is estimated from Eq. (fTT|) . We can then calcu- 
late the number of particles corresponding to the given 
Woods-Saxon profile. If it deviates from the desired 
value by more than SA, we reject the current values 
of po,Pgas,o- and continue. If it is close to the desired 
value within 5A, we calculate also the standard devia- 
tion, ai oop . If this is lower than the stored value of a, 
then a is assigned the new value ai oop and po, p ga s, a are 
assigned the corresponding new values. And so on, un- 
til all values have been tried. It can happen that no 
combination of values satisfies the criterion (|A4j) . The 
reason may be, e.g., an anomalous decrease or increase 
of the gas density close to the cell radius, due to the 
boundary conditions. Then the Hartree-Fock calculation 
is rejected. This explains the seemingly random gaps in 
our microscopic results along isotopic chains. 
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